On the Numerical Stability Analysis of Pipelined Krylov Subspace Methods

Algebraic solvers based on preconditioned Krylov subspace methods are among the most powerful tools for large scale numerical computations in applied mathematics, sciences, technology, as well as in emerging applications in social sciences. The study of mathematical properties of Krylov subspace methods, in both the cases of exact and inexact computations, is a very active area of research and many issues in the analytic theory of Krylov subspace methods remain open. Numerical stability issues have been studied since the formulation of the conjugate gradient method in the middle of the last century, with many remarkable results achieved in the years since. Inexact computations in Krylov subspace methods, either due to floating point roundoff error or intentional action motivated by savings in computing time or energy consumption, have two basic effects, namely, slowing down convergence and limiting attainable accuracy. Although the methodologies for their investigation are different, these phenomena are closely related and cannot be separated from one another. As the name suggests, Krylov subspace methods can be viewed as a sequence of projections onto nested subspaces of increasing dimension. They are therefore by their nature implemented as synchronized recurrences. This is the fundamental obstacle to efficient parallel implementation. Standard approaches to overcoming this obstacle described in the literature involve reducing the number of global synchronization points and increasing parallelism in performing arithmetic operations within individual iterations. One such approach, employed by the so-called pipelined Krylov subspace methods, involves overlapping the global communication needed for computing inner products with local arithmetic computations. Recently, the issues of attainable accuracy and delayed convergence caused by inexact computations became of interest in relation to pipelined Krylov subspace methods. In this contribution we recall the related early results and developments in synchronization-reducing Krylov subspace methods, identify the main factors determining possible numerical instabilities, and outline approaches needed for the analysis and understanding of pipelined Krylov subspace methods. We demonstrate the discussed issues numerically using several algorithmic variants of the conjugate gradient method. The paper concludes with a brief perspective on Krylov subspace methods in the forthcoming exascale era.